A hybrid adaptive large neighborhood search for time-dependent open electric vehicle routing problem with hybrid energy replenishment strategies

As competition intensifies, an increasing number of companies opt to outsource their package distribution operations to professional Third-Party Logistics (3PL) fleets. In response to the growing concern over urban pollution, 3PL fleets have begun to deploy Electric Vehicles (EVs) to perform transportation tasks. This paper aims to address the Time-Dependent Open Electric Vehicle Routing Problem with Hybrid Energy Replenishment Strategies (TDOEVRP-HERS) in the context of urban distribution. The study considers the effect of dynamic urban transport networks on EV energy drain and develops an approach for estimating energy consumption. Meanwhile, the research further empowers 3PL fleets to judiciously oscillate between an array of energy replenishment techniques, encompassing both charging and battery swapping. Based on these insights, a Mixed-Integer Programming (MIP) model with the objective of minimizing total distribution costs incurred by the 3PL fleet is formulated. Given the characteristics of the model, a Hybrid Adaptive Large Neighborhood Search (HALNS) is designed, synergistically integrating the explorative prowess of Ant Colony Optimization (ACO) with the localized search potency of Adaptive Large Neighborhood Search (ALNS). The strategic blend leverages the broad-based solution initiation of ACO as a foundational layer for ALNS’s deeper, nuanced refinements. Numerical experiments on a spectrum of test sets corroborate the efficacy of the HALNS: it proficiently designs vehicular itineraries, trims down EV energy requisites, astutely chooses appropriate energy replenishment avenues, and slashes logistics-related outlays. Therefore, this work not only introduces a new hybrid heuristic technique within the EVRP field, providing high-quality solutions but also accentuates its pivotal role in fostering a sustainable trajectory for urban logistics transportation.


Introduction
Amid the continuous development of the economy and the accelerated pace of urbanization, the demand for urban logistics distribution is growing rapidly [1].Traditional urban logistics predominantly relies on Fuel Vehicles (FV), leading to high fuel consumption, severe noise pollution, and significant carbon emissions, which exacerbate urban environmental degradation [2].In contrast, Electric Vehicles (EVs), characterized by low energy consumption, reduced noise, and minimal pollution, have emerged as an optimal solution for urban distribution and are increasingly recognized as essential tools in transport activities [3].As per International Energy Agency [4], global adoption rates of EVs have tripled in the past three years.Projections from McKinsey & Company [5] suggest that by Year 2030, 10 to 50 percent of new vehicle sales in urban settings will be electric.In alignment with this trend, leading global logistics companies are actively promoting electrified transportation fleets.For instance, FedEx aims to achieve over 50% of its distribution volume using EVs worldwide by 2030 and plans to replace all transport vehicles with electric alternatives by 2040 [6].Similarly, DHL intends to operate 14,000 EVs globally by 2025, representing approximately one-third of its entire fleet [7].These data emphasize the salience of the Electric Vehicle Routing Problem (EVRP) research in shaping urban logistics in response to evolving societal dynamics.
However, EVs face significant limitations in urban distribution compared to traditional FVs.The prominent factors include a limited travel range [8] and long charging durations [9].Urban traffic networks are complex and changeable, with congestion occurring from time to time, leading to varying vehicle speeds that impact the power consumption by the battery and, thus, the vehicle's mileage.Consequently, estimating the remaining real-time range during EV journeys poses a substantial challenge to distribution route planning [10].While charging en route can mitigate range anxiety, the current charging technology extends EV charging time greatly beyond those required for refueling FVs, significantly impacting delivery efficiency.Hence, the battery swapping strategy becomes an effective solution to offset the lengthy charging time.By establishing dedicated battery swapping stations, logistics fleets can replace EV batteries within 10 minutes [11], markedly reducing energy replenishment durations.With the gradual improvement of urban charging infrastructure, transportation fleets have begun to adopt different energy replenishment methods according to actual situations.They opt for charging stations when the recharging demand is low to save on distribution costs, while selecting battery swapping stations when the recharging demand is high to enhance delivery efficiency.These different strategies have varying implications for EV distribution costs and efficiency [12].In summary, when optimizing EV urban logistics distribution routes, planners must consider the time-varying factors of the urban traffic network and choose appropriate energy supplement methods for EVs based on different scenarios.
In addition, in today's fiercely competitive market milieu, companies must strive to reduce operating costs to gain a competitive edge.One effective strategy to achieve this target is to hire professional Third-Party Logistics (3PL) fleets for urban distribution [13].The Open Vehicle Routing Problem (OVRP) exemplifies this cooperative mode.A unique feature of OVRP is that after the completion of delivery tasks, the 3PL fleets do not need to return to the enterprise's depot, which greatly reduces unnecessary empty-load travel and effectively reduces distribution costs [14].When devising a solution for open vehicle distribution, replacing traditional FVs with EVs not only saves energy costs but also reduces carbon emissions, aligning with the expectations of social development.Currently, there are few studies on the Open Electric Vehicle Routing Problem (OEVRP), making it highly significant for research.
Based on the aforementioned description, this research integrates authentic factors, including the dynamic nature of urban traffic networks, diverse energy replenishment strategies, and 3PL fleets delivery, into the Electric Vehicle Routing Problem (EVRP).In this context, a novel Time-Dependent Open Electric Vehicle Routing Problem with Hybrid Energy Replenishment Strategies (TDOEVRP-HERS) is proposed to minimize the distribution cost of 3PL fleets paid by the enterprise.Unlike previous research, this study considers elements such as EV self-weight, real-time package load, travel speed, and travel distance when predicting EV energy consumption.It develops methods for estimating EV power consumption and travel time in the dynamic transportation network.To solve the TDOEVRP-HERS, a hybrid metaheuristic combining Ant Colony Optimization (ACO) and Adaptive Large Neighborhood Search (ALNS) is designed.
The main contributions of this study are summarized as follows: (1) TDOEVRP-HERS is studied in this paper, filling a research gap.(2) A method for estimating EV energy consumption in time-varying traffic networks is devised, aiding logistics fleets to better perform EV urban distribution.(3) The study develops a flexible hybrid energy replenishment strategy that improves delivery efficiency and reduces the distribution cost of logistics fleets by selecting the charging or battery swapping option based on different circumstances.(4) The research pioneers the combined use of ACO and ALNS specifically in the context of EVRPs.Despite each algorithm's individual effectiveness in solving various optimization problems, the exploration of their combined application has remained relatively uncharted.This study leverages the potential synergy between these two algorithms, representing an innovative contribution to the existing body of research.The results demonstrate that the Hybrid Adaptive Large Neighborhood Search (HALNS) offers enhanced problem-solving capabilities, thereby establishing a promising new direction in the field of TDOEVRP.(5) The paper evaluates the performance of the proposed algorithm through extensive numerical experiments and examines the economic and environmental advantages of TDOEVRP-HERS.The experimental findings offer valuable insights and recommendations to governmental policymakers and logistics fleet managers from multiple perspectives, facilitating efforts to promote energy conservation and emission reduction within the logistics sector.
The remainder of this article unfolds as follows: The Section 2 encompasses a review of pertinent literature and elaborates on the innovative aspects of this study.The Section 3 introduces the theoretical methodology for estimating EV power consumption and travel time in the dynamic urban traffic network.The Section 4 formally defines the research question and develops a Mixed-Integer Programming (MIP) model for the problem.The Section 5 details the design of the HALNS.In the Section 6, a sequence of numerical experiments is conducted to validate the efficacy and rationality of the proposed algorithm, accompanied by practical recommendations related to sustainable urban distribution.Finally, the Section 7 encapsulates the conclusion of this study.

Related literature review
As a well-known Combinatorial Optimization Problem (COP), the Vehicle Routing Problem (VRP) aims to determine the optimal set of routes for delivering packages from the depot to customers based on various objective functions [15].Currently, researchers have applied VRPs to distinct transportation scenarios, such as urban waste recycling [16], post-earthquake relief distribution [17], and vaccine emergency response [18].This paper delves into urban distribution through the lens of 3PL fleets.The introduced TDOEVRP-HERS, a VRP variant, integrates time-varying vehicle speeds, varied energy replenishment strategies, and open routing.This section offers a literature review concerning these elements and the metaheuristics employed to address such challenges.focused on TDVRP, which involved studying the real-time speed impact of different weather and congestion levels on vehicle routing [19], identifying ways to evade traffic congestion by analyzing urban traffic's time-varying characteristics [20], and addressing the transport of valuable goods by proposing a secure TDVRP with time windows including pickup and delivery with uncertain demands [21].Researchers then discovered that time-varying vehicle speeds caused various carbon emissions due to different fuel consumption, leading to the development of the Time-Dependent Green Vehicle Routing Problem (TDGVRP).Soysal et al. [22] and C ¸imen et al. [23] established the TDGVRP model, which considered time-varying vehicle speeds, fuel consumption, and carbon emissions.TDGVRP was later focused on cold chain transportation by Guo et al. [24], who proposed a two-stage algorithm that enabled vehicles to wait in place after service to avoid harsh traffic conditions.Recently, with the popularity of EV fleets, the Time-Dependent Electric Vehicle Routing Problem (TDEVRP) is rising gradually.TDEVRP is, in fact, a subset of TDGVRP, which considers time-varying vehicle speeds' impact on EVs' energy consumption rather than FVs' fuel consumption.Lu et al. [25] considered the impact of charging decisions and congestion conditions on the overall delivery process based on traditional constraints in TDEVRP.Bi and Tang [26] examined TDEVRP by incorporating time-varying random traffic conditions and employing the analytical battery model to determine EVs' charging and discharging patterns.They aimed to minimize the overall service time and developed a hybrid rollout algorithm to tackle the problem.Zhang et al. [10] and Keskin et al. [27] approached TDEVRP from different viewpoints.The former incorporated charging during congestion periods into the objective function, while the latter considered the impact of queuing time at charging stations on overall distribution.The literature above provides a theoretical framework and solution methodology for TDEVRP.

Research on the Open Vehicle Routing Problem (OVRP)
OVRP has found practical application in a variety of scenarios, such as newspaper distribution [28], express service [29], campus shuttle [30], and home health care [31].Scholars have examined the OVRP from diverse perspectives considering different constraints.Atefi et al. [32] developed a unique OEVRP model that utilized multiple carriers to deliver packages for a Canadian major producer of cookies and snacks.Azadeh et al. [33] introduced a closed-open hybrid VRP based on the cooperation model of self-operated fleets and 3PL fleets, allowing some vehicles to return to the depot, while the rest of the vehicles do not.Niu et al. [34] discussed the OVRP of outsourcing distribution to a 3PL fleet, and demonstrated the advantages of outsourcing logistics through cost composition analysis.Brandão [35] studied a multi-depot OVRP with the goal of minimizing the travel distance of vehicles.Xia and Fu [36] introduced customer satisfaction rate into the OVRP as a constraint, and constructed an optimization model to minimize the total distribution costs.In addition, many scholars have designed various algorithms for OVRP.Faiz et al. [37] presented two integer linear programming models for OVRP.The first one is an arc-based mixed integer linear programming model solved by a general-purpose solver, while the second one is based on a path-based formulation, which solved by a column generation framework.Chen et al. [38] developed an urgency level-based insertion heuristic to construct an initial solution, and employed a reinforcement learning based variable neighborhood search algorithm to solve the OVRP.Brandão [39] proposed an iterated local search algorithm for OVRP, and solved it with test sets including Solomon, Homberger and Gehring.With the rise of green logistics in recent years, the Green Open Vehicle Routing Problem (GOVRP) has gradually attracted attention.Niu et al. [40] considered vehicle fuel consumption and its environmental impact in GOVRP, conducted a simulation experiment based on the actual road conditions in Beijing, and found that compared to closed routes, open routes reduced total cost by 20%, with fuel consumption costs and carbon emission costs lowered by almost 30%.The above literature provides methodological guidance for OEVRP.However, there is a dearth of literature on OEVRP, indicating that related research has not received sufficient attention.

Research on energy replenishment strategies in EVRP
Efficient energy supplement strategies have always been a prominent research topic in EVRP, aiming to alleviate range anxiety and enhance delivery efficiency in EV distribution.Initially, scholars mainly focused on charging strategies, such as full charging, which has been utilized in the studies of Lin et al. [41], Granada-Echeverri et al. [42], and Kucukoglu et al. [43].Full charging strategy involves completely charging the EV before departing from the charging station, which can effectively mitigate range anxiety.However, this strategy incurs high charging costs and long charging durations, may easily leading to a decrease in customer satisfaction.The strategy of partial charging, in its assessment of charging costs and service timeliness, affords greater flexibility in charging operations.Corte ´s-Murcia et al. [44] and Zhou et al. [45] presented partial charging strategies for EVs when designing EVRP charging schemes.The antecedent enables EVs to be dispatched for multiple distribution tasks, whereas the latter deftly employs the partial charging period to facilitate the delivery of parcels to clients closer to the charging station.In addition, there are different charging technologies in the charging process.Schiffer and Walther [46] and Do ¨nmez et al. [47] proposed a linear charging function when designing the EVRP charging scheme.The former considers a single fast-charging mode at the charging station, while the latter considers a combination of three charging speeds: ordinary charging, fast charging, and super-fast charging.In recent years, battery swapping has emerged as a viable solution for energy replenishment in EV distribution.Jie et al. [48] and Zhou et al. [49] discussed the EVRP considering battery swapping strategy.The former developed a hybrid algorithm to tackle the model, while the latter proposed a multi-objective whale optimization algorithm for the model.Moreover, Raeesi and Zografos [50] devised a synchronized mobile battery swapping strategy for the EVRP, involving the replacement of depleted batteries with fully charged ones at designated times and locations by battery-swapping vans.The literature above can offer theoretical references and methodological support for EVs' energy replenishment strategy in distribution.

Research on metaheuristics of solving VRP
The application of efficient meta-heuristic algorithms to optimize routing problems has been a focus of research in logistics transportation [51].The TDOEVRP-HERS poses unique challenges due to the time-dependent nature of its constraints, alongside the specifics of 3PL fleets, such as non-mandatory return trips, flexible energy replenishment options, and uncertain battery level consumption.Two algorithmic strategies that have shown promise in addressing these complexities are the ACO [52] and ALNS [53].ACO is a population-based metaheuristic that draws inspiration from the foraging behavior of ants [54].Its utilization in complex VRPs has received extensive research attention.Su and Fan [55] used ACO to minimize penalty costs in the Green Vehicle Routing Problem (GVRP), incorporating fuel consumption, carbon emissions, and customer satisfaction constraints.Similarly, Zhang et al. [56] utilized ACO to cut down energy usage in EVRP.Despite its successes, the original ACO has shortcomings, leading researchers to develop the improved ACO.Zhang et al. [57] enhanced traditional ACO with three mutation operators for the multi-objective VRP with flexible time windows, improving local search and global exploration.Li et al. [58] adjusted pheromone update methods of ACO, focusing on optimal solutions for the multi-objective multi-depot GVRP.Jia et al.
[59] created a bi-level ACO for EVRP with capacity constraints.The upper level, unconstrained by electricity, generates routes via ACO, while the lower level uses a heuristic to determine charging plans, using its solutions to update the upper level's pheromone.ALNS is a well-regarded metaheuristic that has successfully addressed numerous COPs.It operates by iteratively exploring a set of candidate solutions in the neighborhood of the current solution and adaptively selecting the most promising ones [60].Scholars have employed ALNS to solve many VRP variants.Sun et al. [61] addressed the time-dependent profitable pickup and delivery problem with time windows by ALNS.Chen et al. [62] discussed a VRP with time windows and delivery robots and developed ALNS to tackle it.Kuhn et al. [63] used the ALNS to integrate order picking and vehicle routing for micro-store supply.ACO is known for generating high-quality initial solutions that respect the problem's constraints.ALNS is particularly effective when dealing with dynamic environments, given its capacity to continually adjust and refine solutions [64].The combination of ACO's global exploration and ALNS's local exploitation may lead to better-quality solutions.The comprehensive search process increases the chance of finding optimal or near-optimal solutions for the TDOEVRP-HERS.However, the hybrid application of these two powerful algorithms in the field of EVRP is currently underresearched, suggesting a need for further exploration.

The research gaps
Table 1 presents a structured categorization of the literature pertinent to TDOEVRP-HERS, drawing from the elements outlined previously.This table systematically delineates distinctions among energy replenishment strategies, routing modes, and solution methodologies within existing studies.Its design facilitates an efficient understanding of the prevailing stateof-the-art and highlights potential research gaps or avenues addressed by this research.
The existing literature lays a solid groundwork for the in-depth investigation of EVRP in urban distribution.However, a review of the literature has identified several gaps in research that merit further exploration: (1) Many studies limit EVs to using charging stations as their sole energy replenishment way during distribution.However, in reality, EV fleets have started to adopt a variety of energy replenishment strategies in current urban distribution.As one of the world's largest consumers of EVs, China explicitly outlined its intention in the "New Energy Vehicle Industry Development Plan (2021-2035)" to "accelerate the construction of battery swap infrastructures," indicating that battery swapping will become a primary method of supplementing EV energy in the future.Despite some scholars have explored EVRP with battery swapping, there remains a significant shortage of EVRP research that considers hybrid energy replenishment strategies with both charging and battery swapping.(2) The time-varying traffic conditions of urban road networks lead to EVs traveling at varying speeds throughout the day, with speed being the most significant factor affecting both EV travel time and power consumption.The fluctuating vehicle speed makes calculating EV energy consumption and remaining mileage challenging, thereby complicating the planning of urban distribution routes.Currently, the academic community has not given sufficient attention to developing reasonable rules for predicting EV energy consumption during urban distribution.(3) To reduce operational costs, an increasing number of enterprises are outsourcing their logistics transportation business.Logistics outsourcing involves a 3PL fleet that departs from the enterprise's depot to perform distribution tasks without having to return to the depot upon completing the tasks, known as open route distribution.Although this cooperation mode is gaining popularity in urban distribution, there is a dearth of literature on TDOEVRP, which warrants further exploration.To fill these research gaps, this study investigates TDOEVR-P-HERS, considering customer demand, service time, EV capacity, time-varying speed, different energy supplement strategies, and open routing.A MIP model is constructed to minimize the total distribution cost of transport fleets, and a HALNS is designed to tackle this model.The study's results provide decision-making references for companies interested in hiring environmental-friendly 3PL fleets and governments seeking to promote sustainable urban logistics.

Calculation of EV power consumption and travel time in a dynamic transportation network
The road network in urban areas is characterized by its time-varying nature, meaning that the speed of vehicles may fluctuate during various time periods when traveling on urban roads.The time taken to traverse the road segment (i,j) depends not only on the road distance but also on the EV's starting time, traveling speed, and real-time load.The vehicle may require several time periods to travel from node i to node j.Therefore, this study refers to the method proposed by Liu et al. [20] to partition a day of 24 hours into multiple equal time periods.H represents the length of a time period, and n represents the number of time periods in a day,  [65] and Basso et al. [8], the energy consumption e R ijk generated by EV k when traveling at speed v R ijk on road segment (i,j) during the R-th time period can be calculated as follows: where ϕ d denotes the output efficiency parameter of the driving motor, φ d denotes the output efficiency parameter of the battery.Other factors influencing EV energy consumption include the rolling resistance coefficient (C r ), air resistance coefficient (R c ), air density (ρ), and the windward area of the vehicle (A w ).θ ij is the slope of road segment (i,j), and g is the gravitational acceleration.L denotes the self-weight of the EV, and u ijk denotes the remaining load of EV k traveling on road (i,j).t R ijk is the travel time of EV k on road section (i,j) during the R-th time period.
Let d ij denote the distance between node i and node j.Let d R ijk denote the remaining distance that EV k still needs to travel on road segment (i,j) at the end of the R-th time period.Under the time-varying urban road network, the following are procedures to calculate the power consumption E ijk and travel time T ijk that EV k requires to traverse road segment (i,j): Procedure 1: The power consumption and travel time generated by EV k during the R-th time period when it leaves node i.Let @ R ik denote the moment when EV k departs from node i within the R-th time period, i.e., T RÀ 1 � @ R ik � T R .Therefore, the feasible travel time of EV k during the R-th time period is T R À @ R ik .Consequently, the distance traveled by EV k on road segment (i,j) during the R-th time period is it implies that EV k has covered the distance between node i and node j within the R-th time period, hence Therefore, the energy consumption of EV k during the R-th time period is e R ijk .Proceed to Procedure 2. Procedure 2: Computation of EV k's power consumption and travel time on road segment (i,j) after the R-th time period.Step 1: Set ξ = 1.Step 2: Calculate the possible travel distance F Rþx ijk of EV k during the (R+ξ)th time period, where , it means that EV k still cannot travel to node j in the (R+ξ)th time period.Consequently, the travel time of EV k during the (R+ξ)th time period is t Rþx ijk ¼ H, and the remaining distance to node j is

Problem and mathematical model
In this section, the TDOEVRP-HERS is formally defined, followed by a presentation of the mathematical formulation that captures the key objectives and constraints of the problem.Factors including energy replenishment strategies, time-dependent vehicle speed, open routing, real-time load, fixed EV usage costs, vehicle travel costs, customer service costs, and energy costs are considered.The goal is to minimize the aggregate of these costs by efficiently dispatching EVs from 3PL fleets and devising optimal vehicle routes.

Problem framework
This subsection introduces and analyzes the TDOEVRP-HERS in urban distribution.A retail company, lacking its own self-operated transportation fleets, delegates all package distribution tasks to a professional 3PL fleet.This 3PL fleet employs homogeneous EVs to deliver identical products to customers scattered all over a city.The objective is to minimize the total operating costs for the 3PL fleet, utilizing up to K vehicles.A complete graph G = {V,A} is defined where V = {0,1,. ..,n} is the set of all vertices, 0 stands for the depot, and A = {(i,j):i,j2V and i6 ¼j} is the set of arcs connecting nodes.Each arc (i,j)2A has an associated distance d ij .Unlike other studies, the infrastructure for EV energy replenishment in this study is bifurcated into two types: Charging stations (F) and battery swapping stations (B).So, Each customer i2C is assigned one and only one EV to fulfill their delivery requirement, with the customer possessing a non-negative demand D i and a designated service time s i .Specifically, before departing from the retail company's depot, all EVs within the 3PL fleet load the maximum load capacity of W kilograms in goods for distribution.Concurrently, the EV fleet is fully charged and primed for action.During the transportation process, priority is accorded to EVs for customer service.However, should an EV find its energy insufficient to sustain the task, it must undertake an energy replenishment.According to the remaining delivery tasks, the EV can choose to replenish the battery level at either a charging station or a battery swapping station.Once the energy replenishment is completed, the EV resumes serving the remaining customers.Upon completion of all transportation tasks, the EV fleet can proceed to other enterprises for additional distribution tasks, without having to return to the retail company's depot.The schematic of the TDOEVRP-HERS is depicted in Fig 1 .In this problem, all the necessary parameters, such as EVs' maximal load W and battery capacity Q, depot coordinates, customers' location, delivery volume, and service time, are known.The study aims to optimize the 3PL fleet's transportation routes, while meeting all customers' distribution demands, and minimizing the total distribution costs that the retail enterprise has to pay.Distribution costs include EV fixed dispatch costs, travel costs, customer service costs, as well as charging and battery swapping costs.Please see S1 Appendix for all the parameters and decision variables used.

Assumptions
To clearly construct the MIP model, the ensuing assumptions are defined: (1) Each EV in the 3PL fleet is dispatched at most once, and the delivery volume of each customer is less than EV capacity.(2) Urban traffic network affects EV speed, which varies with time.(3) The energy expenditure of EVs is affected by a gamut of factors encompassing vehicular self-weight, travel velocity, real-time load, and traversed distance.(4) Once an EV enters a charging station or battery swapping station, it can directly replenish energy without having to wait in line.To ensure the timeliness of distribution, the logistics fleet limits each EV to enter a charging station only once, whereas the number of times it enters a battery swapping station is unlimited.
(5) The charging station adopts a full charging strategy with a constant charging rate, while the battery swapping station has a fixed battery swapping time and the battery level of EVs becomes full upon completion of battery swapping.

TDOEVRP-HERS model
Based on the problem description and variables mentioned earlier, the optimization model of TDOEVRP-HERS is built below: The objective function (2) minimizes the total distribution cost for the 3PL fleet.This cost includes the fixed vehicle dispatching cost (P 1 ), travel cost (P 2 ), customer service cost (P 3 ), and energy cost (P 4 ).Notably, P 4 is the sum of charging and battery swapping costs.Expressions for these costs within the framework of TDOEVRP-HERS are presented in Formulas (3) to (6).In tandem with objective function (2), the study addresses various constraints, organized into five distinct groups: arc and vertex constraints, vehicle capacity constraints, battery level constraints, time constraints, and binary constraints.The arc and vertex constraints are presented through constraints ( 7) to (10), while vehicle capacity constraints are articulated by constraints (11) and (12).The battery level constraints are captured by constraints ( 13) to (18), and time-related restrictions are laid out in constraints (19) to (22).Lastly, constraint (23) enforces a binary constraint.
Constraint (7) stipulates that each EV is dispatched at most once.By enforcing this restriction, the same EV is not redundantly assigned to multiple delivery routes, thereby exerting an indirect effect on the costs associated with vehicle dispatch (P 1 ).Constraint (8) necessitates that each client is visited uniquely once.In the context of the objective function, this constraint aids in minimizing the travel and service costs (P 2 and P 3 ), as it discourages unnecessary detours or repeat visits.Constraint (9) enforces flow conservation among customer locations, battery swapping stations, and charging stations.This constraint indirectly impacts energy costs (P 4 ) by influencing the distribution of EVs and their energy replenishment.Constraint (10) states that each EV is allowed to be recharged at most once on the route.Regarding the objective function, this constraint directly affects the energy cost component (P 4 ) by placing a cap on recharging events, which can influence the allocation of EVs and their corresponding energy consumption patterns.
Constraint (11) asserts that the duty payload apportioned to an individual EV must not transgress its supreme capacity.From the perspective of the objective function, this constraint influences the travel and energy costs (P 2 and P 4 ) by constraining the goods that each EV can carry.Inefficient or overloaded routes that could potentially inflate costs are avoided.Constraint (12) states that the load of an EV after serving a client is determined by subtracting the client's demand from the load before service.This constraint directly impacts the customer service cost (P 3 ), as it influences the load adjustments after each service operation.
Constraint (13) mandates every EV to be completely charged prior to its departure from the depot.This constraint directly affects the energy cost (P 4 ), as it ensures that the initial energy state of each EV is at its maximum capacity, thus minimizing the need for immediate recharging or swapping.Constraint (14) expounds that the EV battery level remains unchanged at customer points.This constraint impacts the energy cost (P 4 ) by maintaining a consistent energy level throughout the service process, which in turn affects the EV's overall energy consumption.Constraint (15) upholds the energy equilibrium of each EV during its transit from the preceding node i to the subsequent node j, considering the energy consumption along the arc (i,j).This constraint directly influences the energy cost (P 4 ) by accounting for the energy consumed during travel and distribution tasks.Constraint (16) ensures that EV's remaining energy is not negative when it reaches any node.This constraint affects the feasibility of the solution and indirectly impacts the EV dispatching cost (P 1 ) by preventing scenarios where an EV's energy level becomes insufficient for completing the assigned tasks.Constraint (17) defines the amount of charging at a charging station, which directly affects the energy cost (P 4 ) by determining the energy replenishment needed during the route.Constraint (18) indicates that when the EV departs from the battery swapping station, its battery level becomes full.This constraint impacts the energy cost (P 4 ) by ensuring that the EV commences its journey with a fully charged battery after swapping.
Constraint (19) embodies the correlation amongst the EV's arrival time, service duration, and departure time at a customer vertex.This constraint affects both the travel cost (P 2 ) and the customer service cost (P 3 ) in the objective function.Constraint (20) signifies the temporal interdependence between the EV's exit time from the present node i and its entry time at the ensuing node j.This constraint impacts the travel cost (P 2 ) within the objective function, as it makes the EV's route go through different traffic conditions.Constraint (21) calculates EV's charging duration at a charging station.This constraint directly affects the energy cost (P 4 ) within the objective function, as it influences the charging-related expenses incurred during the distribution process.Constraint (22) signifies the time relationship before and after the EV replenishes energy at a charging station or battery swapping station.This constraint impacts the objective function by regulating the timing of energy replenishment, thereby affecting both the energy cost (P 4 ) and the overall distribution efficiency.
x ijk ; w ik ; y ik ; z ik 2 f0; 1g Constraint (23), a binary constraint, introduces a discrete decision-making element into the optimization framework.The binary constraint adds a layer of restriction to the optimization process, as the decision to activate or deactivate certain components directly affects the distribution costs captured by the objective function.Depending on the problem specifics, constraint (23) might influence the allocation of vehicles to particular routes, the utilization of charging or battery swapping stations, or the selection of specific distribution tasks.

Algorithm design
The ALNS metaheuristic was originally introduced by Pisinger and Ropke [66] and has since been adapted to various optimization problems.For instance, Koch et al. [67] employed an ALNS framework for an integrated vehicle routing problem with complex loading constraints involving pickups and deliveries.In a similar vein, Z ˇulj et al. [68] devised a hybrid ALNS that integrates tabu search components for tackling an order batching problem.Notably, Seydanlou et al. [69] tailored a multi-neighborhood search algorithm to address sustainable closedloop supply chain management in Iran's agricultural sector.Moreover, Fathollahi-Fard et al. [70] applied an ALNS approach to the generalized quadratic assignment problem and introduced an efficient Benders reformulation based on reformulation linearization technique inequalities.Given this versatility, the ALNS metaheuristic presents a promising foundation for the TDOEVRP-HERS.Expanding upon this foundation, the current framework is extended into a HALNS, drawing inspiration from the ACO metaheuristic.Subsection 5.1 elaborates on the general structure of the HALNS metaheuristic, followed by Subsections 5.2 to 5.5, which provide detailed insights into the distinct constituents of the complete HALNS algorithm.The proposed HALNS ingeniously integrates the ALNS technique with the ACO method, synergizing their respective strengths.The ALNS, with its capability for intensive local search, is coupled with the ACO's global exploration potential.This fusion is achieved by initializing the solution process using ACO's probabilistic approach, generating diverse initial solutions that serve as a foundation for ALNS to perform more refined local searches.The adaptability of ALNS complements the ACO's explorative prowess, enabling it to fine-tune solutions and guide the search towards promising regions of the solution space.This harmonized approach leverages the ALNS to enhance the solutions found by ACO, ultimately contributing to improved convergence and solution quality.The implementation framework of HALNS is meticulously outlined in Algorithm 1.The process begins by designing an ACO to generate a higher-quality initial solution.Next, five destruction operators-Random Removal (RdR), Basic Worst Removal (BWR), Related Removal (RlR), Single Point Removal (SPR), and Station-based Removal (SbR)-are introduced.In addition, three repair operators are included: Improved Random Insertion (IRdI), Improved Greedy Insertion (IGI), and Improved Regret Insertion (IRgI).During the ALNS iteration, these operators are randomly selected to work on the current solution.Simultaneously, an appropriate acceptance criterion is established to adaptively adjust the operators' weights and update the current and optimal solutions.The termination conditions for the algorithm may include the maximum number of iterations or the longest runtime.

The construction of the initial solution
The quality of the initial solution significantly influences the ensuing optimization of the ALNS, and a superior initial solution increases the likelihood of achieving a satisfactory solution.Given the characteristics of the TDOEVRP-HERS, the study employs the ACO to procure an improved initial solution.First, all ants are allocated to depot 0. Each ant then embarks on its journey.Based on the transition probability Formula (24), ants employ a roulette-wheel selection process to choose the customer nodes they visit.These selected customer points are added sequentially to the current route until the ants are no longer eligible to continue customer visits.At this point, the ants are reset to depot 0 and initiate their journey anew.This operation is repeated until all customers in the area have been visited by the ants.
where P m ij is the transition probability for ant m to travel from node i to node j. unvisit m is the set of customer nodes that ant m has yet to visit.visited m is the set of customer nodes that ant m has already visited.τ ij is the pheromone heuristic parameter, and ϑ ij is the expectation heuristic parameter, ϑ ij = 1/d ij .β 1 and β 2 are the weights for pheromone and expectation heuristic factors, respectively.Once all ants have completed visiting all customers, a global pheromone update is executed on the current travel plan, comprising all ants, as per Formula (25).
where t old ij is the pheromone before update, while t new ij is the updated pheromone.rho is the pheromone evaporation rate with 0�rho<1.Set rho = 0.2.Dt m ij signifies the pheromone increment on arc (i,j) for ant m. f is a constant that represents the amount of pheromone secreted by an ant in each travel.Set f = 5.Distance m is the total distance traveled by ant m.
The specific steps of the ACO are as follows: Step 1: Parameter initialization.Begin by inputting the necessary test data, which include the location coordinates of all nodes in the transport network, customer delivery volume, EV's maximum load and battery capacity, congestion time, and time-dependent vehicle speed.Let Iter Ant represent the current iteration, set Iter Ant = 1.Let maxIter Ant represent the maximum number of iterations, set maxIter Ant = 200.Let minLen represent the shortest distance among all ants' travel, set minLen = Inf.Let antNum represent the number of ants, set antNum = 30.Let initialRoute represent the initial optimal travel route.
Step 2: Determine whether the Iter Ant satisfies the relationship of Iter Ant �maxIter Ant .If it satisfies this condition, assign ant m with m = 1, and proceed to Step 3 (Here, tabu m is the ant m's taboo list, and tabu m = ;); otherwise, advance to Step 11.
Step 3: Determine whether m satisfies m�antNum, and if so, enter Step 4; otherwise, skip to Step 9.
Step 4: Place ant m in depot 0 and assign it the maximum allowable load.
Step 5: Based on the transition probability (Formula 24), select a customer node that fulfills the constraints and add it into the current solution.Then, remove this node from ant m's unvisit m set, and place it in ant m's tabu m set.
Step 6: Repeat the operation of Step 5 until ant m no longer meets the requirements to visit any of the remaining customer nodes.
Step 7: Reposition ant m to depot 0 and reset ant m's maximum load.Determine whether the unvisit m is an empty set.If unvisit m 6 ¼ ;, revert to Step 5; otherwise, save ant m's total travel distance antLen m to the iterLen iter set, and save ant m's travel routes antSol m to the iterSol iter set, proceeding subsequently to Step 8.
Step 8: Let m = m+1, clear the taboo list tabu m , reset the unvisit m set, and then proceed to Step 3.
Step 9: Analyze all travel distances within the iterLen iter set to derive the optimal travel distance bestLen iter and the optimal travel routes bestRoute iter .If bestLen iter �minLen, set minLen = bestLen iter and update the globally optimal routes as initialRoute = bestRoute iter .Conversely, if bestLen iter >minLen, then set bestLen iter = minLen and bestRoute iter = initialRoute.Subsequently, pheromones along the route are updated following Formula (25) before progressing to Step 10.
Step 11: Following the principle of greed, insert the charging station or battery swapping station into the optimal initial route, initialRoute, according to the battery level constraints.This process culminates in the acquisition of the initial solution, initialSol, tailored for ALNS.

Operator design
In this study, a comprehensive set of destruction and repair operators are employed within the HALNS.In addition to the fundamental RdR, BWR, and RlR operators, the algorithm introduces two novel operators, SPR and SbR, within the destruction phase, thereby contributing originality to the approach.Moreover, to enhance the algorithm's performance, three improved repair operators-IRdI, IGI, and IRgI-are adapted from the existing literature.This judicious combination facilitates the efficient restoration of solutions.The sequence of operations unfolds as follows: five destruction operators-RdR, BWR, RlR, SPR, and SbRdisassemble the currentSol, which is then reconstructed using the repair operators-IRdI, IGI, and IRgI.Throughout this iterative process, redundant charging or battery swapping stations within currentSol are methodically removed, and new energy replenishment strategies are reintegrated during the repair phase to address route infeasibility.(2) BWR: The cost of a customer node i in the current solution is denoted as (3) RlR: Initially, a customer node i is randomly selected for removal and placed in the C unvisit set.The similarity function is then used to calculate the similarity between each customer node in the current solution and the removed node.The node with the highest similarity in the current solution is selected for removal and added to the C unvisit set.Subsequently, a customer node is randomly selected from C unvisit , and the similarity calculation process is repeated, leading to the removal of the node with the highest similarity.This step is repeated until r customer nodes have been removed.The similarity function used in this paper is defined as follows: where α 1 and α 2 represent weight coefficients satisfying the condition α 1 +α 2 = 1.Set α 1 = 0.6, α 2 = 0.4.d ij signifies the distance between two customer nodes, and |D i −D j | represents the difference in delivery volume between two nodes.If nodes i and j are on the same route, then η ij = 5; otherwise, η ij = 0.
(4) SPR: SPR is an approach that involves the defined areas between two energy replenishment stations, or between a replenishment station and a depot, known as a station service area.This process initiates by randomly selecting a service area and a customer node within that area as an operational point.Subsequently, all customer nodes either to the left or right of the operational point within the service area are deleted.The number of removed nodes is counted as r d , leaving a remaining removal capacity, denoted as r a , where r d �r and r a = r −r d .If r a >0, another station service area is randomly selected and the same removal operation is performed, with a simultaneous update of r a .In the final round of the removal operation, if a situation where r d >r a arises, r a nodes are randomly chosen for removal from the set of r d nodes.Algorithm 5 illustrates the implementation framework of the SPR operator.(5) SbR: SbR begins with the random selection of either a charging station or a battery-swapping station.From the route associated with this selected station, the nearest m customer nodes are removed, where m can be divisible by r.This procedure is repeated until r customer nodes have been removed.The study sets m = 10.Algorithm 6 illustrates the implementation framework of the SbR operator.

Repair operator.
Given that the basic random repair operator, greedy repair operator, and regret repair operator do not account for the impact of customer node re-insertion on route feasibility, this paper proposes three advanced repair operators: IRdI, IGI, and IRgI.Each repair operator involves a two-stage restoration process.The first stage focuses on ensuring the feasibility of customer delivery routes, while the second stage attends to the feasibility of vehicle energy consumption.
(1) IRdI: From the removed customer set C unvisit , randomly select one customer and identify all permissible insertion locations based on load constraints.Subsequently, an insertable position is randomly selected for customer insertion.Repeat this step until all r customers have been reintegrated into the solution.Subsequently, routes that fail to meet power constraints, as determined by time-varying factors and energy consumption estimates, undergo repair via the insertion of energy replenishment stations.Algorithm 7 illustrates the implementation framework of the IRdI operator.(2) IGI: Initiate by randomly selecting one customer from the set of removed customers, C unvi- sit .Identify all feasible insertion positions for this customer, considering load constraints.Subsequently, compute the objective value addition resulting from the insertion of the customer at these potential locations.Choose the location that contributes the least added value for customer insertion.Continue this procedure until all removed r customers have been reintegrated into the solution.Afterward, assess the routes considering time-dependent factors and energy consumption estimates, and incorporate energy-supplementary stations into routes that do not satisfy power constraints for the necessary adjustments.Algorithm 8 illustrates the implementation framework of the IGI operator.
(3) IRgI: Initially, identify all possible insertion locations for the removed r customers, in accordance with load constraints.Next, calculate the additional objective value incurred by inserting customers at different feasible positions.Define customer i's regret value, Re i , as where VC 1 i represents the least objective value added post customer i insertion, and VC 2 i signifies the second smallest additional objective value following the insertion of customer i.Compare the regret values of the removed r customers, and choose the one with the largest regret value for insertion at the location with the least objective value increased.Repeat this process until r customers have been reintroduced into the solution.Ultimately, accommodate any routes failing to meet power constraints by incorporating charging or battery swapping stations, drawing upon time-varying considerations and energy consumption estimates for guidance.

Acceptance criteria
The acceptance criteria for this research utilize the simulated annealing approach as proposed by Adulyasak et al. [72].During the iterative process, if the updated solution updateSol demonstrates better performance than the current solution currentSol, the updated solution is retained.Otherwise, the updated solution is maintained with a probability of p ¼ e À ðzðupdateSolÞÀ zðcurrentSolÞÞ=T .The initial temperature is set to T 0 = 10000 and follows a cooling schedule defined by T n = cT n−1 , where the cooling rate c = 0.995.The entire search process concludes once the number of iterations reaches a predefined maximum limit.

The specific realization of the adaptive process
The fundamental principle of ALNS involves the adaptive adjustment of an operator's weight based on its historical performance.Initially, all operators are assigned equal weight and score, with this study setting the initial operator weight to 10 and the score to 0. Throughout the algorithm's iterative process, the destruction and repair operators are first selected following the roulette wheel rule.Subsequently, these operators are scored differently based on the quality of the updated solution (updateSol) after each iteration.When the number of iterations attains the pre-set threshold, the weights of the operators are updated.As a result, operators with greater weights are more likely to be selected in subsequent iterations.
The search process of ALNS employed in this study consists of N a stages, each of which must execute N b iterations.Thus, the maximum number of iterations, maxIter = N a �N b .Set N a = 4, N b = 50, i.e., maxIter = 200.Assume the weight of operator i in stage j as ω ij and the usage probability of operator i in stage j as p ij .Then, under the current weight, p ij can be calculated by the formula hj , where H represents the set of operators to which operator i belongs.Define the score of operator i in stage j as ε ij .In this study, four levels of scores are established, and the update rules for ε ij are as follows: In the scoring scheme, the corresponding destruction and repair operators accrue points based on the quality of the updated solution relative to both the global optimal and the current solution.Specifically, if the updated solution surpasses the global optimal solution, the operators involved receive an additional 50 points.If the updated solution, while not surpassing the global optimal, improves upon the current solution, the operators gain 20 points.In cases where the updated solution is not an improvement over the current one but is still accepted, the operators are credited with 10 points.Lastly, if the updated solution neither improves upon the current solution nor is accepted, no points are added to the corresponding destruction and repair operators.
The update rules of the operator weights are as follows: where ω i(j+1) denotes the weight of operator i in stage (j+1), θ is the weight adjustment coefficient, indicating the importance of historical weight and operator performance when the operator weight is updated.Set θ = 0.3.Additionally, π ij represents the number of times operator i is selected during the N b iterations of stage j.Should operator i not be utilized during the current stage j, its weight remains unchanged in the subsequent stage.Following a weight update, both ε ij and π ij are reset to zero.

Experiments and analysis
This section offers a detailed account of the test works conducted on the proposed model and algorithm.It covers data preparation, algorithm parameter tuning, experimental evaluations of the superiority of the TDOEVRP-HERS and HALNS, and a thorough analysis of the computational experiences with the HALNS.

Data description
Currently, there is no standardized test dataset specifically tailored for TDOEVRP-HERS, and the geographical distribution of customers varies when 3PL fleets undertake transportation tasks for different enterprises.Consequently, this study refines Solomon's VRP test sets [73], encompassing clustered distribution test sets (C-type), random distribution test sets (R-type), and randomly clustered distribution test sets (RC-type), utilizing these improved test instances as the experimental dataset.Notably, numerous scholars have generated new test instances based on Solomon's test sets.For example, Schneider et al. [74] employed the Solomon test sets as a benchmark and introduced 21 charging stations into each test instance, thereby creating the EVRP experimental dataset they necessitated.Building upon the work of Schneider et al., this study further improves the test instances by converting some of the charging station nodes into battery swapping station nodes.Therefore, new test sets tailored for TDOEVRP-HERS are constructed.Each test instance consists of one depot, 100 customers, five charging stations, and five battery swapping stations.The data provided in these test instances include the coordinates of all vertices, as well as the delivery volumes and service durations for each customer.
To meet the test requirements of this paper, the following data are supplemented: (1) Referring to the update frequency of the traffic congestion index in Beijing, set the duration of each time period as H = 15 minutes, thereby dividing the whole day into 96 time periods.(2) The earliest working time for the enterprise's depot is set as 7:00 a.m.During the traffic peak hours of morning (8:00-9:00) and evening (17:30-19:00), the third-party EV fleet operated at a speed of 20 km/h on the road.(3) The EV fleet is assumed to operate with three time-varying speeds during normal time periods.Utilize the remainder function χ = mod(R,3) to determine these three different velocity values, where R represents the R-th time period.When χ takes the value 1, 2, or 0, respectively, it corresponds to the respective time-varying travel speeds of the EV as 54 km/h, 72 km/h, or 42 km/h.( 4 The prescribed methodologies within the HALNS framework are implemented using the MATLAB R2020b software and performed on a microcomputer, furnished with a 3.60 GHz processor and a memory capacity of 16 GB of RAM.

Parameter tuning
An initial series of experiments is undertaken to fine-tune the parameters of the ALNS metaheuristic, as detailed in Table 2.For this purpose, 27 representative instances, as proposed by Ticha et al. [75], are chosen.These include the SOL instances r101, r105, c103, c104, rc101, and rc105, each with 50 customers, and the NEWLET instances 1, 2, and 3, each comprising 50 customers and 100 nodes.These nine instances are considered for three correlation tiers: No-Correlation (NC), Weak Correlation (WC), and Strong Correlation (SC).
Utilizing these instances, the following procedure is embraced.Initially, the parameters of the RlR heuristic are tuned.The ALNS scheme is then implemented, limited to this destruction operator and the IGI operator.Successively, attention is directed to one of the parameters, testing several values for it.For each value, the tuning instances are solved five times; the value that consistently manifests the best average solution quality is ultimately adopted.
For the SbR heuristic, an identical methodology is applied.Other parameters are systematically adjusted in a similar manner, but the comprehensive ALNS scheme is employed rather than relying on individual destruction and repair heuristics.

Experimental evaluation of the proposed model and algorithm
This section aims to validate the superiority of the proposed model and algorithm through five distinct experiments.Each experiment corresponds to a separate subsection and focuses on a specific aspect of the model or algorithm.Two experiments, found in Subsections 6.3.1 and 6.3.5, are dedicated to examining the effectiveness of the HALNS algorithm.Subsection 6.3.1 evaluates the efficiency of the HALNS in solving test sets featuring diverse geographic distributions of customers.Subsection 6.3.5 assesses the HALNS's robustness to environmental disturbances, with a specific emphasis on its ability to achieve high-quality solutions amidst varying levels of traffic congestion.The other three experiments, detailed in Subsections 6.3.2, 6.3.3, and 6.3.4,are designed to validate the advantages of the TDOEVRP-HERS model.Subsection 6.3.2 investigates whether switching from FVs to EVs indeed leads to a reduction in the overall transportation costs of logistics fleets while also promoting a sustainable urban environment.Subsection 6.3.3 analyzes whether a hybrid energy replenishment strategy outperforms a single replenishment strategy.Finally, Subsection 6.3.4 examines whether employing a 3PL fleet is more economical and environmentally friendly than utilizing a self-operated fleet.In this way, the five experiments comprehensively evaluate the performance of the HALNS algorithm and the TDOEVRP-HERS model, providing valuable insights for practitioners and policymakers.

EV route planning for different test instances.
Various test instances, featuring diverse customer geographic distributions, serve to ascertain the feasibility of the HALNS proposed herein.Table 3 shows the test results, where TN signifies the test set name, TC signifies the total distribution cost incurred by the 3PL fleet (in yuan), FC signifies the fixed EV dispatching cost (in yuan), SC signifies the customer service cost (in yuan), DC signifies the vehicle travel cost (in yuan), EC signifies the EV energy replenishment cost (in yuan), EN signifies the number of energy replenishment times for all EVs, DN signifies the number of EVs enabled for distribution, RT signifies the running time of HALNS (in second), and AVE denotes the average value.
The test results depicted in Table 3 reveal that: (1) The fixed EV dispatching costs, travel costs, and customer service costs account for an average of 29.27%, 39.48%, and 18.30% of the total distribution costs, respectively, totaling 87.04%.This finding suggests that vehicle fixed dispatching costs, travel costs, and customer service costs remain the primary factors affecting distribution costs in urban distribution.In practice, when enterprises engage 3PL fleets for urban distribution, they should require the fleet to use EVs with larger capacities to minimize the number of EVs dispatched.At the same time, to avoid increased travel costs due to traffic congestion, logistics fleets should consider the effect of time-varying vehicle speeds on travel time when performing distribution tasks.Furthermore, it is recommended that fleet managers train fleet members to reduce customer service time and improve service efficiency, which can not only reduce the service cost that enterprises have to pay for the fleet, but also improve the competitiveness of the transport fleet.These efforts can ultimately reduce the total distribution cost of the fleet and foster stronger collaborations between enterprises and 3PL fleets.(2) EV energy replenishment costs during distribution account for a range of 6.91% to 19.20% of the total distribution cost, with an average of 12.82%.It is noteworthy that the majority of EVs on distribution routes only require a single energy supplement at most to complete transportation tasks.Thus, using EVs for urban distribution has a negligible effect on distribution efficiency.Compared to traditional FVs, EVs incur lower energy consumption and offer additional advantages such as zero emissions and reduced noise levels.Embracing EVs for urban distribution can aid logistics fleets in reducing energy consumption and gaining cost edges, all while decreasing the environmental impact of transport activities.This aligns with the vision of policymakers who seek to promote green and sustainable development of urban distribution.While FVs do not account for energy supplement time during distribution, they generate fuel consumption and carbon emissions.Following the methodology and experimental parameters employed by Liu et al. [20] for calculating fuel consumption and carbon emissions of FVs in their VRP study, the associated costs are incorporated into the distribution expenses of the logistics fleet.Table 4 presents the results of comparative experiments between TDOEVRP and TDOFVRP (Time-Dependent Open Fuel Vehicle Routing Problem), where TD signifies the total EV travel distance, GC signifies fuel consumption costs, and CC signifies carbon emission costs incurred by FVs.Remaining symbols maintain their previously established definitions.
The test results depicted in Table 4 reveal that: (1) Using EVs instead of FVs for urban distribution offers significant cost savings for logistics fleets.Specifically, substituting EVs for FVs leads to an average reduction of 31.03% in distribution costs, with the energy consumption cost for the logistics fleet decreasing from an average of 36.56% to 12.36% of the total distribution cost-a remarkable decline of 66.20%.This finding highlights the advantages of EV fleets in urban distribution, including energy savings and reduced operational expenses.To enhance competitiveness, logistics fleets should adopt EVs for their distribution operations, which aligns with society's call for sustainable development within the logistics industry.Government departments should promote EV distribution by strengthening the construction and management of EV energy replenishment infrastructures and enhancing the convenience and safety of EV fleet use.(2) Using EVs for urban distribution leads to a mere 4.69% average increase in total distance traveled compared to using FVs.This finding can be attributed to the well-developed energy replenishment infrastructure in urban areas, which makes it easy for EVs to locate nearby charging stations or battery swapping stations.Additionally, using FVs for urban distribution generates an average carbon emission cost of 14.60 yuan, while EVs produce no carbon emissions, improving urban air quality and promoting sustainable logistics transportation.However, it is worth noting that the carbon emission cost of FV distribution accounts for only 0.92% of the total distribution cost, which is less than 1% of the proportion, indicating that logistics fleets do not prioritize reducing emissions when using FVs for urban distribution.This result suggests that China's current carbon price is too low to incentivize logistics fleets to curtail emissions actively.Therefore, Chinese policymakers should formulate more reasonable carbon trading-related policies to promote sustainable logistics transportation.

Comparative tests of different energy replenishment strategies.
A comparative experiment is conducted to assess the effects of various energy replenishment strategies on EV distribution costs, including Hybrid Energy Replenishment Strategy (HERS), Pure Charging Strategy (PCS), and Pure Battery Swapping Strategy (PBSS).The experiment uses multiple types of test sets while keeping the other parameters unchanged.The comparative test results of the three energy replenishment strategies are presented in Table 5, where TT is the total distribution time for all EVs (in minutes), which encompasses EVs' travel time on the route, customer service time, and energy replenishment time, and ET is the total energy replenishment time for all EVs (in minutes).Remaining symbols maintain their previously established definitions.
The test results depicted in Table 5 reveal that: (1) In each test instance, the EV energy replenishment time achieved with PBSS is significantly lower than the other two strategies.Compared to HERS, PBSS reduces EV energy supplement time by a range of 51.30% to 69.25%, with an average saving of 58.20%.Compared to PCS, PBSS exhibits even more substantial savings in EV energy supplement time, with an impressive average time-saving of 79.80%.However, it should be noted that PBSS comes with higher distribution costs for logistics fleets.The total distribution costs incurred by PBSS are the highest among all test instances.On average, the use of PBSS increases the distribution cost of the logistics fleet by 4.74% compared to PCS.Therefore, when the 3PL fleet undertakes distribution tasks for enterprises, it is crucial to assess the urgency of customer package delivery.If there are time-sensitive customers, opting for battery swapping stations for EV energy replenishment can enhance delivery efficiency.On the other hand, if customers do not prioritize timeliness, charging stations can be utilized to reduce distribution costs.(2) Although HERS does not offer the optimal solution in terms of both total distribution cost and total distribution time, it exhibits a modest increase in distribution cost of only 2.65% on average compared to PCS, which has the lowest distribution cost.Similarly, compared to PBSS, which achieves the shortest distribution time, HERS shows a slight increase in distribution time by an average of 2.94%.Despite PCS being cost-effective, it requires an average of 5.69% more distribution time than HERS.Likewise, although PBSS reduces distribution time, it incurs an average of 1.90% higher distribution costs than HERS.These findings highlight the ability of HERS to strike a balance between total distribution costs and distribution efficiency, providing flexibility in energy replenishment operations during transportation tasks.Therefore, using HERS to support the urban distribution of EV fleets is a more scientific and reasonable approach.

Comparative tests of open routing and closed routing.
A comparative experiment is conducted to analyze the differences in urban distribution route planning between self-operated transport fleets and 3PL fleets.Multi-type test sets are utilized while keeping the relevant parameters unchanged.It should be noted that the self-operated transport fleet follows a closed routing, which requires returning to the enterprise's depot upon completing the distribution task.Conversely, the 3PL fleet adopts an open routing, eliminating the need to return to the depot after task completion.The test results are presented in Table 6, where TCSR denotes the proportion of total distribution cost saved by hiring a 3PL fleet instead of a self-operated transport fleet (in %), and TDSR denotes the proportion of total distribution distance saved by hiring a 3PL fleet (in %).Remaining symbols maintain their previously established definitions.The test results depicted in Table 6 reveal that: (1) Engaging a 3PL fleet leads to an average savings of 6.20% and 5.61% in total distribution cost and time, respectively, compared to deploying a self-operated transport fleet.This advantage stems from the fact that when a thirdparty EV fleet is hired, EVs are not required to return to the company's depot after completing the distribution task.As a result, the empty travel time of EVs is reduced, leading to decreased travel costs in logistics transportation.Therefore, for companies with an imperfect warehouse network layout, specifically lacking multiple depots, opting for a third-party EV fleet instead of a self-operated fleet can effectively reduce their operational expenses.(2) Hiring a 3PL fleet leads to a considerable reduction in EV travel distance compared to dispatching a self-operated transport fleet, with an average decrease of 9.35%.This finding highlights the potential for further optimizing EV distribution-energy supplement route planning when employing 3PL fleets, resulting in more efficient distribution solutions for EVs.Therefore, promoting this collaboration model is recommended for practical logistics distribution.6.3.5.Robustness analysis of HALNS under different traffic conditions.While keeping other experimental parameters constant, five distinct combinations of traffic conditions (Ⅰ, Ⅱ, Ⅲ, Ⅳ, Ⅴ) are designed.These conditions are outlined in detail in Table 7, where STC represents severe traffic congestion during the current time period, MTC represents mild traffic congestion during the current time period, TFS indicates that the traffic flows smoothly during the current time period, and the values provided represent the average travel speeds of EVs under the respective traffic scenarios (in km/h).Employing the RC107 test set for experimentation, the test results are presented in Table 8, where DT denotes the travel time on the routes (in min), and PU denotes the total energy consumption of all transport EVs (in kWh).Remaining symbols maintain their previously established definitions.
The test results depicted in Table 8 reveal that: (1) The enhancement of traffic conditions yields a gradual reduction in the travel time of EVs during transit.Transitioning from traffic scenario Ⅰ to scenario Ⅴ, the en-route travel time of EVs diminishes by 42.62%.This finding Ⅲ.This is attributed to the fact that the transition from transportation scenario Ⅲ to scenario Ⅴ results in a 75.87% increase in energy replenishment costs, whereas the reduction in travel costs amounts to only 16.81%.The rate of increase in energy costs surpasses the decrease in travel costs.Thus, it becomes apparent that when planning the transportation route for the urban distribution fleet, fleet managers must establish a reasonable travel speed range to strike the optimal balance between efficiency and total cost.(3) In comparison to traffic scenario Ⅲ, which boasts the lowest total distribution cost, the average disparity between the total distribution cost of other traffic scenarios and scenario Ⅲ is 3.46%.This minor gap in distribution costs showcases the adaptability and route optimization capabilities of the HALNS based on real-time traffic conditions and time-dependent constraints.When traffic conditions deteriorate, the proposed HALNS prioritizes the delivery of customers in closest proximity to minimize travel time.Conversely, when traffic conditions improve, the HALNS adeptly schedules the insertion of energy-supplementary stations, effectively striking a balance between energy costs and congestion impacts.

Algorithm comparative test
To assess the efficacy of the proposed HALNS, a series of comparative analyses are performed, spanning two experimental categories: comparisons with an exact solver and other metaheuristics.For small-scale computations, outcomes from the HALNS are juxtaposed with those from the widely-recognized CPLEX solver, while for large-scale computations, HALNS's optimized results are set against those from two renowned metaheuristics.These benchmarks facilitate a thorough appraisal of the HALNS's performance across varying scales, underscoring its aptitude and preeminence in addressing large-scale optimization challenges.6.4.1.Comparison with exact solver.In order to validate the proposed mathematical model and assess the performance of the HALNS algorithm, the CPLEX solver version 12.6 is used to address small-scale test instances.These results are then juxtaposed against those produced by the HALNS.To ensure the generality of the experiment, selected subsets from Goeke's dataset [76] are utilized, with certain charging station coordinates adjusted to batteryswapping station coordinates to accommodate experimental requirements.Each test instance is run through HALNS 10 times to discern the optimal result.The findings are tabulated in Table 9.Within the table, NC denotes the number of customers, NS indicates the number of charging stations, NB signifies the number of battery swapping stations, and 'gap' highlights the gap between the ALNS-derived solution and that from CPLEX.Remaining symbols maintain their previously established definitions.Table 9 delineates that the HALNS is adept at efficiently tackling small-scale TDOEVRPs-HERS.Given the lower number of vertices, HALNS consistently identifies exact solutions that coincide with those derived from the CPLEX solver.However, considering the NP-hard nature of TDOEVRP-HERS, the computational duration for CPLEX escalates dramatically with increasing instance size.As a result, CPLEX's average execution time for all test instances clocks in at 105.83 seconds, in stark contrast to HALNS's mere average of 2.26 seconds.This observation accentuates HALNS's proficiency in navigating towards optimal solutions while preserving exceptional efficiency.

Comparison with other metaheuristics.
To further validate the superiority of HALNS, an Ant Colony System (ACS) and a Genetic Algorithm (GA) are developed to tackle the TDOEVRP-HERS.The parameter settings of ACS align with the work studied by Li et al. [77], while the parameter settings of GA follow the experimental setup of Park et al. [78].Multi-type test sets are used to compare the performance of HALNS, ACS, and GA.Table 10 reports the test results, with the symbols retaining their previously defined meanings.
The test results depicted in Table 10 reveal that: (1) In each test instance, HALNS consistently outperforms ACS and GA in terms of total distribution cost and travel distance.On average, HALNS achieves savings of 3.37% in distribution cost and 3.18% in travel distance compared to ACS, and 5.09% in distribution cost and 4.31% in travel distance against GA, respectively.In addition, the energy replenishment stations selected in the HALNS solution are conveniently located near the served customers, while the solutions generated by ACS and GA allocate some energy supplement stations far away from customers to the transport fleet, suggesting that ACS and GA may not reasonably optimize the energy supplement route in transportation, which will affect distribution timeliness and lower customer satisfaction.These findings exhibit the effectiveness of HALNS in optimizing logistics distribution routes, improving logistics distribution efficiency, and reducing logistics distribution costs.(2) Despite the strong optimization capabilities of the HALNS, it demands an average of 17.58% more running time than the ACS and 15.52% more than the GA.Fig 6 highlights the time-consuming calculation burden of the three algorithms through a line graph when processing different test instances.This result can be attributed to three key aspects of the HALNS: operational complexity, depth of search, and adaptivity overhead.Firstly, HALNS integrates multiple destroy and repair heuristics, enabling it to explore a broader solution space and thereby increasing its potential to uncover high-quality solutions.However, the utilization of these operators tends to be computationally taxing, and thus extends the overall processing time.Secondly, the adaptive nature of the HALNS means that it will often perform a more thorough and exhaustive search of the solution space, as compared to ACS or GA.This deeper exploration can yield better results but requires more processing time.Finally, unlike simpler algorithms like ACS or GA, HALNS incorporates an adaptivity mechanism where it modifies the use of different operators based on their performance.This feature enhances the flexibility and robustness of HALNS, but it contributes to the increased running time as the algorithm needs to constantly evaluate and adjust the weights assigned to various heuristics.Consequently, optimizing the solving efficiency of HALNS, without compromising solution quality, will significantly influence its broader application in future contexts.

Sensitivity analysis
Upon validating the superior performance of the proposed HALNS, sensitivity analyses are performed on the destruction and repair operators to assess their efficiency and influence on computational costs.For these evaluations, test set R101 is specifically chosen, and the average changes between the initial and final solutions, as determined by the HALNS algorithm, are calculated.These findings are tabulated in Table 11.To scrutinize the efficacy of each operator, the HALNS is executed using a singular pair of destruction-repair operators, and improvements in the solutions are documented.A subsequent analysis, focusing on the individual performance of the destruction and repair operators, processes the average changes noted in Table 11 for each operator and visually represents them in Fig 7.
Results from Table 11 indicate that the optimal pair of heuristics is the combination of RlR with IGI and IRgI.These two pairings yield an average enhancement of approximately 120% from the initial to the final solution.Upon examining  In conclusion, the operators employed in the proposed HALNS are pivotal in shaping computational costs.They enhance search ability by adeptly selecting and integrating elements during the solution formulation.By carefully choosing the most fitting destruction and repair operators, the HALNS can traverse the solution space with greater efficacy and attain marked enhancements in solution quality with diminished computational demands.

Discussion
Optimization algorithms have evolved remarkably, fostering advancements across various domains, ranging from online learning, scheduling, and transportation to medicine and data classification.Given the increasing complexity of contemporary decision problems, incorporating advanced algorithms will be instrumental in deriving more effective solutions.The decision problem tackled in this research, the TDOEVRP-HERS, predominantly relied on the HALNS algorithm.While HALNS offers significant merits, it becomes imperative to juxtapose it with other advanced optimization methods to foster comprehensive solution frameworks.
Advanced optimization techniques, such as hybrid heuristics, metaheuristics, hyperheuristics, self-adaptive algorithms, island algorithms, and polyploid algorithms, bring forth  specialized capabilities, enabling researchers to address the nuanced facets of complex decision problems.Their adaptability and customizability render them suitable for a plethora of applications.For instance, an online-learning-based evolutionary many-objective algorithm has been manifested to optimize decision-making in dynamic environments [79].Adaptive polyploid memetic algorithms have demonstrated effectiveness in scheduling applications, such as optimizing truck schedules at cross-docking terminals [80].In the realm of VRP, algorithms that accommodate both exact and metaheuristic methodologies have been applied, ensuring robust solutions even in multi-objective settings [81].
Given the myriad applications of these algorithms, their relevance to TDOEVRP-HERS cannot be understated.Hybrid meta-heuristic algorithms, for instance, have shown considerable promise in supply chain network designs under uncertain conditions [82], suggesting potential adaptability for electric vehicle routing under time-varying traffic networks.Similarly, the customizability offered by hyper-heuristics can be of paramount importance, especially when devising metaheuristics for continuous optimization [83].This points towards their potential utility in fine-tuning the solutions for TDOEVRP-HERS.
In conclusion, while the HALNS has proven its merit in this research, the vast landscape of advanced optimization algorithms holds significant promise.By embracing a multifaceted approach that incorporates insights from various algorithms and domains, the research community can pave the way for more affluent, nuanced solutions to intricate decision problems like TDOEVRP-HERS.

Conclusions
The study investigates the TDOEVRP-HERS as applicable to urban distribution.Therefore, an approach is devised to estimate EV energy consumption within a dynamic urban traffic network.In focusing on third-party logistics distribution, a MIP model is developed with the objective of minimizing the 3PL fleet's total distribution costs.Based on the model's characteristics, a HALNS is designed to tackle the model.Comprehensive numerical experiments are conducted to verify the feasibility of the proposed model and algorithm.The key findings of the experiments are as follows: 1. Cost components: Principal determinants affecting transportation fleets' distribution costs encompass fixed vehicle dispatching, travel expenses, and customer service costs.Employing EVs in urban transit settings does not lead to exorbitant costs.In fact, EVs are more energy-efficient than FVs, enhancing the competitive edge of 3PL fleets.Their low-carbon footprint further accentuates urban logistical transport's sustainability, mitigating greenhouse gas emissions and conserving the urban environment.
2. 3PL fleets' flexibility: The capacity of 3PL fleets to provide adaptable distribution solutions negates the necessity to return to the depot post-task, optimizing route planning and reducing travel distances.Leveraging the benefits of energy conservation and low emissions by EVs, they can substantially diminish greenhouse gas emissions, endorsing eco-friendly transportation practices.Especially given the rising inclination towards low-emission zones in cities, collaboration with third-party EV fleets becomes essential for enterprises.
3. Energy supplement strategies: The battery swapping method, while enhancing distribution efficacy, entails higher costs.Conversely, the EV charging approach is cost-effective but time-intensive.The hybrid strategy, though not optimal in terms of cost and time, exhibits minimal cost variance when compared to the charging method and mirrors the batteryswapping technique in distribution time.Recognizing the diverse energy demands of varying routes and schedules, hybrid strategies grant the 3PL fleet the leeway to customize energy top-ups.
Based on the above findings, recommendations emerge.Enterprises, especially those grappling with non-optimized warehousing networks, should actively contemplate synergies with third-party EV fleets.Not only does this present an avenue for operational cost reductions, but it also ushers in a shift towards eco-friendly urban logistics.
While this research offers a robust foundation, it is circumscribed to the dynamics of urban distribution, focusing primarily on 3PL distribution scenarios.Thus, broader applications may require additional nuances.Future directions will investigate the potential applicability and modifications of the HALNS algorithm in other logistical scenarios and delve deeper into how destruction-and-repair operators in various urban settings might improve the route planning of 3PL fleets.

Algorithm 2 . 2 :Algorithm 3 .
For a visual depiction, refer to Fig 2, where the leftmost image presents currentSol, the central image demonstrates the outcome of applying the destruction operator (resulting in the removal of six customers and three energy replenishment stations), and the rightmost image showcases updateSol following the repair operator's intervention.5.3.1.Destruction operator.The destruction operator selectively removes a predetermined number of customer nodes from the current solution.Given a total customer count, n, the quantity r designated for removal is determined by the equation r = [α×n].In this equation, α represents the removal ratio, and the [] is the rounding symbol.Set α = 0.1.The customers removed through this operation are subsequently added to the unvisited customer set, C unvisit , while any invalid routes and superfluous charging or battery swapping stations are discarded.The general structure of destruction operators is offered in Algorithm 2. The general framework of destruction operators input: A currentSol and the number of requests to be removed r output: A partial_solution 1: Initialize removal list (C unvisit ;) Use roulette rules to randomly choose a destroy heuristic based on their weights 3: while termination criteria are not met do 4: Apply the selected destruction operator to remove a request r 5: C unvisit C unvisit [{r} 6: end while (1) RdR: Randomly remove r customer nodes from the current solution.Algorithm 3 illustrates the implementation framework of the RdR operator.The general framework of the RdR operator

( 3 ) 6 . 3 . 2 .
The HALNS's running time ranged from a minimum of 211.47 seconds to a maximum of 281.32 seconds, with an average of 255.23 seconds.This result indicates that the HALNS proposed in this study can provide high-quality EV distribution and energy replenishment route planning solutions that meet decision objectives in a reasonable period, demonstrating both high efficiency and feasibility.Fig 3 illustrates the EV distribution and energy replenishment route planning solutions for the C103, C203, R103, and RC203 test instances.The solutions for each test instance involving 100 customers are shown clearly and distinctly in the figure, with few instances of route detours and intersections.Based on test cases reflecting a realistic delivery scale, the proposed algorithm demonstrates its potential to consider various practical factors and offer valuable guidance for the transportation route optimization of logistics fleets.Comparative tests of TDOEVRP and TDOFVRP.A comparative experiment is conducted to analyze the differences in urban distribution costs between the utilization of EVs and FVs.Multi-type test sets are employed while keeping relevant parameters unchanged.

Fig 3 .
Fig 3. EV route planning solutions for different test instances.https://doi.org/10.1371/journal.pone.0291473.g003 Fig 4 visually illustrates the gaps in the objective function values among the three algorithms, clearly showing HALNS routinely yielding the most favorable objective function values.Meanwhile, Fig 5 showcases the EV distribution route schemes generated by HALNS, ACS, and GA for test set C208.Unmistakably, the HALNS-generated solution delineates the most coherent distribution routes with minimal situations of route detours and intersections.

Fig 7 ,
it is evident that RlR serves as the most effective destruction operator (Fig 7(A)), while IGI stands out as the superior repair operator (Fig 7(B)), showcasing a 12.20% overall improvement when compared to IRdI.

Table 1 . A summary of relevant works.
T 1 �; ½T 1 ; T 2 �; :::; ½T nÀ 1 ; T n �g, where [0,T 1 ] represents the first time period of a day, and [T n−1 ,T n ] represents the n-th time period.Additionally, [T R−1 ,T R ] signifies the R-th time period of a day, where R2{1,2,. ..,n−1,n}, and T R −T R−1 = H.The traveling speed of EV k on road segment (i,j) during the R-th time period is denoted as v R ijk .By combining the EV energy consumption calculation methods proposed by Goeke and Schneider The energy consumption of EV k during this time period is e Rþx ijk .Let ξ = ξ+1 and proceed to Step 2. If F Rþx ijk � d RþxÀ 1 ijk, it means that EV k can travel to node j during the (R+ξ)th time period, then the travel time of EV k during this time period is t Rþx ijk ¼ d RþxÀ 1 ijk =v Rþx ijk , and the energy consumption of EV k during this time period is e Rþx ijk .Proceed to Procedure 3. Procedure 3: Computation of the total power consumption and travel time for EV k traversing road segment (i,j).E ijk ¼ X Rþx x¼R e x ijk , T ijk ¼ X Rþx x¼R t x ijk .The calculation process is complete.
Randomly choose a route path D to destroy, and randomly remove a customer node o from path D .
[71], f(s) represents the objective value of the current solution, and f −i (s) is the objective value after removing customer node i.Initially, BWR calculates the objective value of the current solution.Subsequently, it computes the objective values after the removal of different customer nodes.After that, it derives the cost associated with each customer node.Finally, it selects the customer node with the highest cost for removal and repeats above step until r customer nodes have been removed.This strategy was proposed by Hemmelmayr et al.[71].Algorithm 4 illustrates the implementation framework of the BWR operator.

Table 7 . Average speeds of EVs in different traffic conditions.
While the amelioration of traffic conditions boosts the transportation efficiency of 3PL fleets, it also leads to a notable escalation in energy consumption.Shifting from traffic scenario Ⅰ to scenario Ⅴ, the energy consumption of EVs rises by 190.26%.Heightened energy consumption necessitates more frequent energy supplements.Consequently, EVs fail to achieve the lowest total distribution cost in traffic transport scenario Ⅴ, but instead achieve it in transportation scenario